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ABSTRACT 

Observations of Keplerian disks of masers in NCG 4258 and other Seyfert galaxies can be used to 
obtain geometric distance estimates and derive the Hubble constant. The ultimate precision of such 
measurements could be limited by uncertainties in the disk geometry. Using a time-dependent linear 
theory model, we study the evolution of a thin initially eccentric disk under conditions appropriate 
to sub-pc scales in Active Galactic Nuclei. The evolution is controlled by a combination of differ- 
ential precession driven by the disk potential and propagating eccentricity waves that are damped 
by viscosity. A simple estimate yields a circularization timescale of Tdrc ~ 10^(r/0.1 pc)^/^ yr. Nu- 
merical solutions for the eccentricity evolution confirm that damping commences on this timescale, 
but show that the subsequent decay rate of the eccentricity depends upon the uncertain strength of 
viscous damping of eccentricity. If eccentricity waves are important further decay of the eccentricity 
can be slow, with full circularization requiring up to 50 Myr for disks at radii of 0.1 pc to 0.2 pc. 
Observationally, this implies that it is plausible that enough time has elapsed for the eccentricity of 
masing disks to have been substantially damped, but that it may not be justified to assume vanishing 
eccentricity. We predict that during the damping phase the pericenter of the eccentric orbits describes 
a moderately tightly wound spiral with radius. 

Subject headings: accretion, accretion disks — masers — galaxies: active — galaxies: Seyfert — 
galaxies: distances and redshifts — galaxies: individual (NGC 4258) 



1. INTRODUCTION 

Observations of water maser emission from a geomet- 
rically thin disk wi thin the central pc of NGC 4258 
(|Mivoshi et al.|[l995l ) have permitted a geometric deter- 
mina t ion of the distance to th e galaxy (|Herrnstein et al.) 
Il999t iHumphrevs et al.l l2008f ) , and an estimate of the 
Hubble constant. The recent discovery of apparently 
similar maser emission in a number of more distant galax- 
ies has raised the exciting prospect of a percent level de- 
termination of Hq independent of both distance ladder 
uncertainties and other cosmological parameters. Such a 
determination would have considerable cosmological util- 
ity, perhaps most importantly by improving the ability 
of current and future datasets to cons train the behavior 
of dark energy (IHull2005l: IOlling|[200l . 

The precision with which Hq can be determined from 
measurements of masers in a single system is limited by 
uncertainties in the galactic peculiar velocity, observa- 
tional errors, and any unmodeled complexity in the ge- 
ometry of the source. The first two of these limits present 
instrumental and observational challenges, and consider- 
able progress is guaranteed given sensitive, multi-epoch 
observations of maser disks at higher redshift than the 
one in NGC 4258. The third requires both observations 
and theory. Typically, it is assumed that the maser emis- 
sion originates from a thin, circular disk, and this sim- 
plest assum ption is consistent with all existing data on 
NGC 4258 (IHumphrevs et al.l[200l . There is no ob- 
vious theoretical rationale, however, for assuming that 
the eccentricity is exactly zero, or for extrapolating the 
NGC 4258 results to other systems. Indeed, the obser- 
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vation that the masing disk in NGC 4258 is measur- 
ably warped suggests the possibility that different disks 
might exhibit a range of geometries. Although the ge- 
ometry can be constrained observationally, this requires 
the introduction of additional parameters whose fitting 
inevitably degrades the precision with which the distance 
can be measured. 

The purpose of this paper is to outline, from a the- 
oretical perspective, the expected evolution of the ec- 
centricity of a geometrically thin accretion disk at sub- 
pc distances from a supermassive black hole. Using a 
time-dependent linear theory model of an eccentric disk 
(Goodchild & Ogilvic 2006) I examine both the rate at 
which eccentricity decays, and the transient structure of 
the eccentric disk before it becomes circular. 

2. DESCRIPTION OF THE DISK EVOLUTION 

Consider a thin, flat disk in which the vertical scale 
height h <ti r. We define a complex eccentricity via. 



E = e exp [iw] 



(1) 



where e is the magnitude of the eccentricity and w the 
longitude of pericenter of the fluid streamlines. We as- 
sume that at the radii of interest the potential is domi- 
nated by the Keplerian potential of the black hole <f>BH, 
to which must be added an additional contribution due 
to the disk, 

$ = $BH + *disk- (2) 

The fluid in the disk has density p, pressure p, angular 
velocity and adiabatic exponent 7. We assume that 
turbulent stresses within the disk act to damp eccen- 
tricity, and parameterize the efficiency of the damping 
via a bulk viscosity ae written in terms of a Shakura- 
Sunyaev (1973) a parameter. Linearization of the two- 
dimensional fluid equations (in the inviscid limit, except 
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for the aforementioned bulk visco sity) yields an evolu- 
tion equation for the eccentricity (jGoodchild &: Ogilvid 

Eool), 
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The first two terms on the right hand side describe pre- 
cession driven by any non-Keplerian part of the disk po- 
tential and by pressure gradients within the disk. These 
terms do not change the magnitude of the eccentricity. 
The third term has the form of a Schrodinger equation 
- it describes waves of eccentricity that propagate radi- 
ally through the disk and are damped by the action of 
viscosity. 

Theoretically, it is expected that accretion disks in Ac- 
tive Galactic Nuclei (AGN) ought to be self-gravitating 
at the relatively large radii where masing occurs, and 
this permits some simplification of the evolution equa- 
tion. We write the surface density as a power-law in 
radius, 

S - Sor'^ (4) 

and compute the disk potential <i>disk from the enclosed 
disk mass 



M, 



disk 



2'!TrYjdr 



(5) 



assuming a spherically symmetric mass distribution. We 
assume that over the relatively narrow radial range for 
which masing occurs the sound speed Cg can be taken 
(approximately) to be constant, and note that p ~ E/2/i 
and that h — Cs/^. The rate of precession due to pres- 
sure gradients is then oc (a constant at the radii where 
masing occurs), while the rate of precession due to the 
non-Keplerian potential is oc Eg, which rises with the 
disk mass"^. Thei: 



3ir ratio. 
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therefore depends upon the disk mass - in low mass disks 
precession will be dominated by the pressure term while 
in higher mass disks gravity will dominate. The latter 
limit is probably appropriate for masing disks. At the 
radii - of the order of a tenth of a pc - where the masing 
is observed it is likely that disks around supermassive 
black holes are self-gra vitating (Kolychalov fc SunvaevI 
T980; Clarke 1988; Shl osman. Begelman fc FrankI llQQOl : 
Goodman. i2003.) . Self-gravitating disks develop spiral 



structure, whose existence in the NCG 4258 disk is 
hinted at by observations of clustering and as ymmetry in 
the m aser emission regions tMaoz 1995; Mao z fc McKei 
[Tpl). For a self-gravitating disk MBn/Mdisk ~ {h/r)~^, 
^ h/r, and pressure effects are negligible. 



\EJ/\E, 



^ These terms have different signs, so it is possible to choose a 
surface density and sound speed profile for which they cancel, o r for 
which the disk precesses as a rigid body. Indeed, IStatlerl (120011 ) de- 
rived a family of potentially long-lived eccentric disk models which 
owe their stability to rigid body precession. We do not consider 
special solutions of this typo here, since there is no reason to expect 
that masing disks would have a structure that leads to vanishing 
differential precession. 



Writing p — pc^ the evolution equation Q then simpli- 
fies to, 
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We study the eccentricity evolution implied by this equa- 
tion in §3. 

2.1. Limitations of this description 

Equation ([3]) is an extremely simple representation of 
the evolution of an eccentric fluid disk. The most obvious 
limitation is the use of a linear equ ation ra t her th an the 
nonlinear formalism developed bv I Ogilvid ()2001[ ). The 
linear equation ought to provide an approximate descrip- 
tion of the dynamics at late times, when the eccentricity 
is small, but will evidently fail if the initial eccentricity 
and / or its variation with radius is large. A second and 
even more important limitation concerns the uncertain 
nature of the eccentricity damping. We have parameter- 
ized the damping via a bulk viscosity - rather than the 
usual shear viscosity used in the Shakura-Sunyaev theory 
of accretion disks - because disks dominated by a Navier- 
Stokes shear viscosity are unstable to growth of eccentric 
modes (Kato 1983; Ogilvie 2001). This is probably not 
a physical problem, rather it reflects the fact that an- 
gular momentum tr ansport mediated by th e magnetoro- 
tational instability (jBalbus fc Hawlevlll998f ). self-gravity, 
or other physical mechanisms cannot be well described 
via a Navier-Stokes shear viscosity. For our purposes we 
simply assume that the stress in the disk acts such as to 
damp eccentricity, parameterize that damping efficiency 
(arbitrarily) via a bulk viscosity, and treat ae as a free 
parameter. 

3. ECCENTRICITY EVOLUTION 

3.1. Analytic estimates 

The terms in equations ([3]) and ^ describing preces- 
sion do not alter the magnitude of the disk eccentric- 
ity. The wave-like term proportional to 7 does change 
the local eccentricity, but preserves a global invariant 
(iGoodchild fc Ogilvid[200l . 

£^= I -r^pn\E\^dr. (8) 



Damping of £^ occurs due to the viscous term at a rate 

The presence of eccentricity waves means that the ra- 
dial profile of E needs to be considered along with the 
temporal evolution. This requires a numerical solution 
of equation ([7]), which we defer to §3.2. For an estimate, 
however, we can define a local damping timescale, 

pn\E\^ 

which is simply the ratio of the right hand sides of equa- 
tions ([5]) and li E varies on the same scale as the az- 
imuthally averaged disk properties then to order of mag- 
nitude \dE/dr\'^ ~ E^/r^ and we find that. 
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This is just the usual "viscous time" of the disk if 
Qfe, the viscosity responsible for damping eccentricity, is 
comparable to the shear viscosity as driving mass in- 
flow^. It is well known that, at radii of 0.1 pc and be- 
yond, the viscous time in AGN disks is extremely long 
osman. Begelman fc Fran"3 IT99Q ). For the specific 
case of NGC 4258, observations of the thickness of the 
masing disk imply a midplane temperature of T ~ 600 K 
and a sound speed ~ 1.5 kms~^ ()Argon et al.l 120071 ). 
consistent with the conditions needed for H2O maser 
emission. Adopting these values, the viscous time is 

-2 
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yr. 



(12) 

^4xl07yry VO-1 Pc/ 
Unless sub-pc disks in AGN remain unperturbed across 
timescales of the order of a Gyr, circularization due to 
viscous evolution will not occur. 

Much faster circularization occurs as an indirect con- 
sequence of differential precession, which for the disks of 
interest is of the order of {h/r)~^ ^ lO'^ times faster than 
viscous evolution. The precession timescale Tp = 2t:/uj, 
where 
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As noted above, we expect disks at radii of the order of 
0.1 pc around AGN to be self-gravitating, in which case 
(|Toomrdll96l 

where E is the disk surface density. Using this relation 
to fix the surface density in terms of the sound speed 
(assumed constant) and disk radius, we find, 

-1 



Tp = 8.2 X lO^Q 



yr, (15) 



,1.5kms"V VO-lpc, 
independent of the black hole mass. This relatively rapid 
precession has an important consequence. Although it 
does not alter e, the differential precession leads to a 
rapid variation of w with radius. This reduces the char- 
acteristic lengthscale over with E varies in equation (fTO|) , 
and allows viscosity to damp the eccentricity much more 
quickly. 

The origin of the gas in the masing disks of Seyfert 
galaxies is not known. One possibility is that the 
disk is established following the infall of relat ively small 
amou nts of low angular momentum gas (Kin g~fc Pringld 
|2007( ). In this scenario, one expects the disk to be ini- 
tially non-circular, with the fluid streamlines describing 
nested aligned ellipses. Partially motivated by such a 
picture, we consider as initial conditions a disk in which 
e 7^ and w = const. After time t, the pericenters 
of ellipses separated by radial distance Ar have angular 
separation, 

doj . , 
Aw = —Art. 16) 
dr 

^ This argument leads to the rule of thumb that an eccentric disk 
ought to circularize "on the viscous timescale" . This is potentially 
misleading, first because there is no direct physical relation between 
cte and as , and, second, because differential precession can lead to 
a situation in which \dE/dr\'^ 3> jr"^ in much less than a viscous 
time, allowing for faster eccentricity evolution. 



Defining the winding scale as the radial separation for 
which Aw = 27r, 



47rQ 



(17) 



This scale decreases inversely with and as a result the 
timescale on which viscosity can damp the eccentricity 
(equation [T0|) scales with the system age as. 



'^damp 
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(18) 



The disk will circularize when t ~ Tdamp- Defining this 
time as the circularization time Tcirc, we obtain. 



r,i,c = [!^^QYI\GM^^YI^a-'l'c;^l'r^l\ (19) 

which displays only a weak dependence on the black hole 
mass and on the (very uncertain) viscosity damping the 
eccentricity. Since masing action is only possible over 
a fairly limited range of disk temperatures the only sig- 
nificant dependence is the almost linear scaling of the 
circularization timescale with radius. Adopting repre- 
sentative values of the parameters the timescale is. 



W = 8.5 X 10^Q2/3 
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(20) 



This simple calculation provides an estimate of when 
enough differential precession will have accumulated to 
start eccentricity damping. If only differential preces- 
sion and damping are considered then subsequent evolu- 
tion will simply wrap up the eccentricity vector into an 
ever-tighter spiral, hastening the decay of f ^ even fur- 
ther. Indeed, very rapid damping due to this process 
was seen in numerical solutions for precessing warped 
disks in the regime where the warp is transmitted vis- 
cously (jArmitage fc Natarajanl 11993 ). In the eccentric 
case, however, inspection of equation ([7]), shows that 
once damping is important (or even beforehand), the 
wavelike term proportional to 7 is also likely to be sig- 
nificant. This term can offset differential precession, and 
affect the rate at which further decay of the eccentricity 
occurs. As we show subsequently, it turns out that while 
equation ([20]) suffices to predict the onset of eccentricity 
damping fairly accurately, the rate of subsequent decay 
varies with much more strongly than would be sus- 
pected based on the rather weak scaling derived above. 

3.2. Numerical results 

To study the eccentricity evolution in more detail we 
have computed numerical solutions to equation ([7]). We 
discretize the equation on a uniform spatial grid be- 
tween Till and Tout, and advance the system in time us- 
ing a second-order scheme that, in the absence of damp- 
ing, preserves the invariant £^ to sufficient accuracy for 
our purposes. We assume that the inner disk is circu- 
lar, and set E — at the inner boundary Tin. The 
correct boundary conditions to apply at Tout are less 
clear, since the outer edge of the masing region may 
reflect either a physical disk edge (defined, for exam- 
ple, by the onset of disk fragmentation) or merely the 
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Fig. 1. — The evolution of the eccentricity at r = 0.15 pc (upper 
panel) and eccentricity invariant {€'^)^^'^ (lower panel) computed 
from numerical solutions to equation jTj. The units on the vertical 
axes are arbitrary. Three runs are shown, all of which assume 
a Q = 1 disk around a a black hole of mass 4 X 10^ f^Q- The 
disk has a constant sound speed of 1.5 km/s and 7 = 1.4. The 
initial eccentricity is a gaussian centered on 0.15 pc with a width 
of 0.05 pc and a constant ro. The solid curves show results for a disk 
extending from 0.05 pc to 0.4 pc, with = 0.1. The dashed curves 
depict the evolution of an identical disk within larger boundaries 
that extend from 0.025 pc to 0.8 pc. The vertical dashed line 
shows the analytic estimate for the circularization time for these 
parameters. The dot-dashed curves show the effect of increasing 
the damping coefficient to ae = 1- 



outer extent of the zone that prod uces maser emission 
(|Neufeld. Malonev fc Conge^ [l99l . Typically we use 
zero-gradient conditions on the components of E (i.e. 



dE^/dr = dEy/dr = at 



We have run models 



with different boundary conditions, and with different 
choices for rout , to gauge the effect of the boundaries on 
the results. 

Figure [T] shows the decay of an initially gaussian ec- 
centricity perturbation, 



e(r) = eo exp 



[r - rof 



(21) 



centered on tq = 0.15 pc with width Ar = 0.05 pc. As 
this is a linear calculation, the value of eg is entirely 
arbitrary. We assume a Q = 1 disk with a constant 
sound speed Cg = 1.5 km/s, which results in a E cx r^'^/^ 
surface density profile. The black hole mass is set to the 
NGC 4258 value of 4 x 10^ M©. The Figure shows results 
obtained from two models with a,, = 0.1 (which differ 
only in the radial extent of the disk), and one model 
with stronger eccentricity damping ae = 1. We plot 
both the eccentricity of the disk at rp, and the integrated 
eccentricity invariant ^f£^ in a form that is proportional 
to the magnitude of the eccentricity. 

The numerical results confirm that the analytic es- 
timate of the damping timescale (equation 1^0)) is rea- 
sonably accurate. For the two runs with a,, = 0.1 we 
find that \/£^ drops by a factor of two after 10 Myr of 
evolution, as predicted. For these runs, however, subse- 



quent damping is relatively slow. It takes approximately 
50 Myr for \/£^ to decay from its initial value by an 
order of magnitude. The corresponding curve showing 
the local value of e displays clearly the wave-like nature 
of the evolution. Although e displays a general decaying 
trend the evolution is oscillatory rather than a monotonic 
decline, with the period of the oscillations varying with 
the radial extent of the disk. Only for strong damping 
(ae ~ 1) is the wave-like behavior suppressed. In the 
strong damping limit the initial decay of the eccentricity 
is not that much faster than in the weakly damped cases 
but, once damping has set in, the subsequent decline of 
the eccentricity to negligible values is rapid. 

We have also studied the effect of different initial con- 
ditions. In general, we find that the decay rate of \/£^ 
is quite robust against changes in the initial conditions. 
The evolution of the value of e at a particular radius, 
on the other hand, can exhibit large transients at rela- 
tively early times if the initial conditions contain large 
gradients in E. Physically it seems likely that an eccen- 
tric disk would settle down to a smooth distribution of E 
shortly after formation, so the gaussian initial conditions 
we have plotted ought to be representative of the likely 
evolution. 

Figure [2] plots the spiral traced by w{r) at different 
times in the weakly damped run (ae — 0.1) with the 
smaller radial extent of the disk (0.05 pc < r < 0.4 pc). 
The initial evolution of the eccentricity vector is driven 
by differential precession, with damping setting in once 
the spiral has become sufficiently tightly wound - for this 
run a value of rdvo/Ar ~ 27r, which is achieved after a 
few Myr, suffices. During the damping phase there is no 
further increase in the winding due to differential preces- 
sion, which is offset by the wave-like term in the evolu- 
tion equation. As a consequence, final damping of the 
eccentricity in the low ae runs occurs relatively slowly. 
Indeed, as shown in Figure [H for the low ae runs the 
sense of the spiral pattern reverses during the damping 
phase roughly periodically, though only the first of these 
reversals is likely to have any physical significance since 
at later times the actual magnitude of the eccentricity 
is quite small. As with the eccentricity oscillations de- 
scribed above, the reversals are largely eliminated in the 
more strongly damped run. 

Based on the above results, it is clear that while the 
eccentricity is decaying the structure of the disk can 
be complex. It is not even possible to predict securely 
whether the eccentricity vector will describe a leading 
or a trailing spiral pattern. At almost all epochs a de- 
scription of the disk in terms of aligned nested ellipses 
is a very poor approximation. A better description is to 
assume a linear radial variation in tu, 



ro(r) = w{ro) -f 



dnj 
dr 



[r - To). 



(22) 



'"0 



As shown in Figure 2, this description is usually reliable 
over at least a factor of two in radius. Typical values of 
the twist in the eccentricity vector, rdvj/dr are in the 
range of tt to 2tt. We believe that fitting a model of this 
kind to observational data would represent a substantial 
improvement over a model with vj = const, though at 
the cost of introducing an additional free parameter in 
the winding angle of the eccentricity vector. 
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t = 105 yr t = 6.3 xlO^ yr t = 2.55 xlO' yr 

Fig. 2. — The radial variation of the angle of pericenter vj of the fluid streamlines, plotted at different times. The run corresponds to the 
solid curves in Fig. [T] The eccentricity vector is wound up by differential precession until damping sets in. As waves subsequently propagate 
through the disk the spiral pattern periodically reverses, as shown in the rightmost panel. The heavy dashed curve in the rightmost panel 
shows a linear fit to ro(r), which is reasonably accurate at most epochs over a limited range (a factor of two) in radius. 



10 - 




-10 - 



10= 10" 10' 

t / yr 

Fig. 3. — The evolution of the twist rdro/dr in the eccentricity 
vector, evaluated at r = 0.15 pc for the standard run (solid curve, 
with parameters corresponding to the same curve in Fig.[T} and the 
run with Oe = 1. The weakly damped case displays approximately 
periodic reversals between leading and trailing spiral patterns in 
the eccentricity vector. These occur on a timescale of the order 
of 10^ yr. The oscillations are largely suppressed in the run with 
stronger eccentricity damping. 



4. DISCUSSION 

Theory provides scant guidance as to the formation 
mechanism and hfetime of sub-pc disks in AGN. The 
long viscous timescale at these radii ~ of the order of a 
Gyr - means that it is unlikely that such disks formed 
from outward viscous expansion of gas that was initially 
much closer to the black hole, but this still leaves sev- 
eral possibilities for feeding gas directly to sub-pc radii. 
One idea is that clouds on low angular momentum or- 
bits are sporadically shredded by the tidal field of the 
black hole, with the debris cooling and settling into a 



disk (|King fc Pringld l2007f ) . Depending upon the disk 
mass and radius involved, the result could be a stable 
gaseous disk or gravitational fragmentat ion into a disk 
of stars ()Gammiell200ll: iRice et al.|[2003[) . Such an idea 
provides one possible explanation for the disk of young 
star s observed at sub-pc radii fro m our own Galactic Cen- 
ter ([Levin fc Beloborodovl 120031 ). If this is the origin of 
the masing disks then it is probable that the initial con- 
ditions would involve eccentric and / or warped disks. 
However, there are also plausible scenarios that would 
yield initially circular disks. For M31, for example, it 
has been suggested that the innermost, circular, disk of 
stars formed from a gravitationally unstable disk that ac- 
cum ulated from the stellar winds of stars situated further 
out (jChang et a l. 2007) . Given these very different mod- 
els, the observation that the mas ing disk in NGC 4258 is 
at least approximately circular (|Humphrevs et al.]|2008f) 
raises two questions. First, does the lack of obvious ec- 
centricity constrain the formation mechanism to favor 
the scenarios that yield initially circular orbits? Second, 
it is reasonable to assume that the (unmeasured) eccen- 
tricity is actually vanishingly small? In this paper we 
have attempted to provide partial answers to both of 
these questions. 

Our main result is that the initial damping of a thin 
self-gravitating disk at sub-pc distances from a black hole 
is relatively rapid. The timescale at 0.1 pc is around 
10 Myr, and even this may be an overestimate since any 
nonlinear effects (such as the formation of shocks as ellip- 
tical orbits precess) will hasten the decay further. Rapid 
damping is driven by differential precession, which occurs 
on a much shorter timescale than the viscous evolution 
of the disk. The lack of gross eccentricity in NGC 4258 
is therefore not surprising, and does not constrain the 
disk formation mechanism. Once damping is underway, 
however, it can take a surprisingly long period for the ec- 
centricity to decay to negligible levels. If the viscous coef- 
ficient governing eccentricity damping is relatively small 
(cte = 0.1 in our description) the presence of waves in 
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the disk allows significant eccentricity to survive for up 
to 50 Myr. Evidently it would be useful to determine the 
efficiency of eccentricity damping under the conditions 
likely to prevail at sub-pc radii, perhaps via improved 
simulations of self-gravitating eccentric disks, which to 
date have focused on simp ler questions such as w hether 
the disk fragments or not ([Alexander et al.|[2008[ ). Cur- 
rently, though, it seems unwise to conclude that a good 
upper limit on the eccentricity implies strictly circular 
orbits. 

Finally, we addressed the question of the predicted or- 
bital structure of an eccentric disk during the circulariza- 
tion process. Existin g constraints on the ecc entricity of 
the NGC 4258 disk (jHumphrevs et al.l l2008l) hav e been 
derived assuming a model, developed bv IStatleH ()200lD 
for different purposes, in which the eccentric orbits are 



aligned and nested. In the case of very thin masing disks 
this model is not a good approximation to the likely 
structure, which instead involves a moderately tightly 
wound spiral in the angle of periapse. A reasonable fit, 
over a modest range of radii, is possible using the next or- 
der approximation in which the periapse angle is a linear 
function of orbital radius. 
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grant AST 0407040. 
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